Clearing the Workspace

rm(list=ls())

Bash Code to Calculate PCA from Extracted SNP’s

#!/bin/bash

activate gill

plink --vcf /dss/dsshome1/09/ra78pec/Local_PCA/chrom5.recode.vcf \
      --keep-allele-order \
      --allow-no-sex  \
      --const-fid \
      --set-missing-var-ids @:# \
      --pca \
      --out /dss/dsshome1/09/ra78pec/Local_PCA/chr5_pca

Plotting Local PCA

From Chromosome Position 35925000 - 37875000

library(dplyr)
library(tidyr)
library(plotly)


axis = list(showline=FALSE,
            zeroline=FALSE,
            gridcolor='#ffff',
            ticklen=4,
            titlefont=list(size=13))


eigenvec <- read.table("/Users/vicegill/Documents/Project_2022/chr5_pca.eigenvec")
eigenvec <- eigenvec[,2:22]
eigenvec <- separate(eigenvec,V2,c("Location","ID"),sep="_")
eigenvec$ID <- as.factor(eigenvec$ID)
fig <- eigenvec %>% plot_ly() %>% add_trace(type='splom',
                                            dimensions=list(
                                              list(label='PC1', values=~V3),
                                              list(label='PC2', values=~V4),
                                              list(label='PC3', values=~V5),
                                              list(label='PC4', values=~V6)
                                            ),color= ~ID,colors = c('#636EFA','#EF553B','#00CC96',"purple"),
                                            marker=list(size=7,line = list(
                                              width = 1,
                                              color = 'rgb(230,230,230,230)'))
                                            )
fig <-  fig %>% style(diagonal = list(visible = FALSE)) %>% layout(
  hovermode='closest',
  dragmode= 'select',
  plot_bgcolor='rgba(240,240,240, 0.95)',
  xaxis=list(domain=NULL, showline=F, zeroline=F, gridcolor='#ffff', ticklen=4),
  yaxis=list(domain=NULL, showline=F, zeroline=F, gridcolor='#ffff', ticklen=4),
  xaxis2=axis,
  xaxis3=axis,
  xaxis4=axis,
  yaxis2=axis,
  yaxis3=axis,
  yaxis4=axis)
fig